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ABSTRACT 

The GMRT Epoch of Reionization (EoR) experiment is an ongoing effort to measure 
the power spectrum from neutral hydrogen at high redshift. We have previously re- 
ported an upper limit of (70 mK)^ at wavenumbers oi k ^ 0.65 /iMpc~^ using a basic 
piecewise-linear foreground subtraction. In this paper we explore the use of a singular 
value decomposition to remove foregrounds with fewer assumptions about the fore- 
ground structure. Using this method we also quantify, for the first time, the signal loss 
due to the foreground filter and present new power spectra adjusted for this loss, pro- 
viding a revised measurement of a 2cr upper limit at (248 mK)^ for k = 0.50 /iMpc~^. 
While this revised limit is larger than previously reported, we believe it to be more 
robust and still represents the most up-to-date constraints on reionization at z ~ 8.6. 
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1 INTRODUCTION 

The Epoch of Reionization (EoR) began as the first stars 
ionized the neutral hydrogen around them, and ended 
when that ionization extended a cross most of the Hub- 
ble sphere. IFurlanetto et al.l (|2006l ) provides a thorough re- 
view of the subject. Based on the electron column den- 
sity to the CMB and under the assumption that reion- 
ization was instantaneous, WMAP data s uggests it would 
have occurred at z = 10.4 (jKomatsu et al.ll2011. V Theoreti- 
cal work, however, has generally sug gested that reionization 
was a patchy and extended process (Furlanetto et al. 2004; 
McQ uinn et al. 2007; Zahn et al. 2007; Friedrich et al. 2011; 
Su et al.ll201ll : lGriffen et al.ll2012l ). Observations of absorp- 



tion lines in quasar spectra can be used to limit the frac- 
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tion of neutral hydrogen at high redshift (jGunn PetersonI 
1 19651 1 and it is generally accepted that reionization was 
complete by a red shift of z 6 (jDiorgovski et al. I I2OOII : 
iBecker et al.ll200ll l, t hough the actual HI fraction may still 
have been quite high (jMcGreer et ahllioill : ISchroeder et al.l 
J012). Using the glo bal 21 cm signal as a function of redshift, 
[Bowman k Roger^ (|201QI 1 have put a lower li mit on t he du - 
ration of the EoR of > 0.06, while .Zahn et al.1 (|2012l ) 
have used measurements of the kinetic Sunyaev-Zel'dovich 
effect with the South Pole Telescope to suggest an upper 
limit on the transition from a neutral fraction of 0.99 to 
0.20 of Az < 4.4. 

The redshifted 21 cm Hi spectral line can be used to 
trace the patchy distribution of neutral hydrogen in the 
universe before t he first luminous source s formed until the 
end of the EoR (Furla netto et al.1 12OO6Y The distribution 
during the transi tion can be used to constrain cosniolog- 
ical parameters (|Furlanetto et al.1 l2009l : iMao et al.1 I2OO8I : 
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Coorav et al] l2008l : iMcQuinn et"all l20Q6l : IPandolfi et aP 
201lh and deduce the nature of the first ionizing sources 



themselves (llliev et al.ll2012l: Kovetz Kamionkowskill2012l : 



iDatta et al.1 I2QI2I : iMaiumdar et al.1 I2OI2I ) , including pos- 
sible exotic reioniza tion scenarios (e.g., IFurlanetto et al.l 
I2OO6I : iHaimanl [20TTI ). Unfortunately the 21cm signal is 
many orders of magnitude less than foregrounds from 
gala ctic and extragalactic sources at the relevant f requen- 
cies (|0h Mackll2003l : lde Oliveira-Costa et al.l l2008^ These 
foregrounds are currently one of the largest obstacles to de- 
tecting the 21 cm signal and several schemes have be en de- 
veloped to address th e problem fe.g., Petrovic Ohi[201ll: 
Liu TegmarkI l2Qlll: IChapman et al.1 l2Q12l : IPillon et al.1 



2OI2I : IParsons et al.ll2012h 



Several groups are ma king progress towards measuring 
the 21 cm power spectrum. Bowma n et al.l (|2008l ) estimated 
an upper limit to the contribution of HI to the redshifted 
21 cm brightness temperature of 450 mK. The Precision Ar- 
ray for Probing the Epoch of Reionization (PAPER) has 
repor ted a limit of appro ximately 5K with a 310 mK noise 
level ('Parsons et al. 2010") . The Murchinson Widefield Array 
(MWA; Lonsdale et al. 2009) is expected to be able to de- 
tect both the amplitude and slope of the power spectrum 
with a signal-to- noise ratio > 10 (|Beardslev et al.l l2012h . 
Both PAPER and MWA have emphasized the importance 
of antenna l ayout in max imizing se nsitivity to t he EoR sig- 
nal (Parsons et al.ll2012l : iBeardslev et al.ll2012l ). The Low- 
Frequency Array ( ^LOFAR) names the EoR as one of its Key 
Science Projects ([Harker et al.ll20l"ol : iBrentiens et aDl2011 ) 
and is c urrentl y bei ng commissioned in the Netherlands. 
IZaroubi et al.1 (|2012l l estimate that LOFAR wih have the 
potential to overcome the low signal-to-noise to directly im- 
age the neutral hydrogen. Future generations of telescopes, 
in particular the Square Kilometer Array (SKA), should be 
capable of direct imaging (Carilli et al. 2004) but will no t 
be in full operation for another decade (|Rawlings Schilizzil 

[mil). 



The GMRT-EoR experiment has been an ongoing ef- 
fort using the Giant Metrewave Radio Tel escope in India 
(|Swarup et al.|[l99ll : lAnanthakrishnanlll995l l, which in con- 
trast to other experiments features large steerable antennas 
with a collecting area com parable to LOFAR , and a rela- 
tively small field of view. In lPaciga et al.l (|201lh we reported 
an upper limit on the neutral hydrogen power spectrum of 
(70mK)2 at 2a using a simple piecewise-linear foreground 
filter. However, this limit did not account for any 21 cm sig- 
nal lost in the foreground filter itself. The purpose of the 
current work is to quantify the potential signal loss, and to 
present these results with a new SVD foreground filter. 

This paper is organized as follows. In section [2] we briefly 
describe the data and the preliminary data analysis. In sec- 
tion [3] we discuss the SVD foreground filter, followed by 
quantifying the signal loss it causes. Finally, we make an es- 
timate of the full 3D power spectrum in section [4l and con- 
clude in section [5l When necessary we will use the WMAP-7 
maximum likelihood paramet ers Qm — 0.271, Oa — 0.729, 
and Ho = 70.3 km s"^ Mpc~^ (iKomatsu et al.ll201lh . All dis- 
tances are in comoving units. 




10^ 

Mode number 



Figure 1. An example singular value spectrum for the three 
shortest baselines, each approximately 50-52 wavelengths at the 
zenith, with the largest singular value normalized to 1 for each. 



2 OBSERVATIONS AND DATA 

The data analysed in this paper were taken over 5 nights in 
December 2007 and total about 40 hours. The observations 
were centred on PSR B0823+26, a pulsar about 30 degrees 
off the galactic plane, which is used to calibrate both the 
phases and the ionosphere by dividing the pulsar period into 
several gates and separating the on- and off-pulse gates. The 
primary beam has a FWHM of 3.1 degrees and a maximum 
angular resolution of about 20 arcsec. The bandwidth cov- 
ers a frequency range from 139.3 MHz to 156.0 MHz in 64 
frequency bins of 0.25 MHz each with a time resolution of 
64 seconds. This corresponds to a redshift range of z = 8.1- 
9.2. For more details on the obser v ations , including RFI re- 
moval strategies, see 'Paciga et alj (|201lh . The remainder of 
this section outlines the differences in the data analysis com- 
pared to this earlier work. 

In addition to automated flagging of visibilities and the 
SVD RFI removal pipeline for broadband interference, we 
have also added manual flagging of faulty antennas, times- 
tamps, and frequency ranges that are exceptionally noisy or 
that visually appear to have RFI left after the automated 
procedures. Approximately 15 per cent of the visibilities are 
flagged in this way, while the dynamic range, defined as the 
ratio of the peak flux of the calibration source to the RMS 
outside the primary beam, is improved by as much as a fac- 
tor of four. 

To improve the comparability of each of the 5 nights 
of observing, we limit each night to the same LST range, 
and have regridded the visibilities in time such that each 
night shares the exact same timestamps, equally spaced in 
one minute intervals. Additionally, it is known that the flux 
of a pulsar can change significantly with time, which creates 
another source of variability from night to night. Since the 
visibilities are gated on the period of the pulsar, it can easily 
be subtracted from the data. These two changes reduce the 
RMS noise in the difference between pairs of days by a factor 
of 2 on average. 

Finally, the power spectrum is calculated from the cross- 
correlation of pairs of nights in annuli of (u, v) space. This 
gives the 2D power perpendicular to the line of sight as 
a function of baseline length |u| or equivalently multipole 
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Figure 2. Sky image before and after an SVD foreground subtraction, for the night of December 10, 2007, using baselines up to 4 km. 
The colour scales are in units of Janskys and the black circle shows the full width half maximum (FWHM) of the primary beam. The 
RMS before any foreground removal (left) is 50mJy. After removing 8 SVD modes (right) the peak goes from 1.6 Jy to 39mJy with an 
RMS of 2 mJy. Residual point sources can still be seen around the edge of the beam while those within about one degree of the centre 
are effectively removed. 
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Figure 3. Sky images of the top 8 SVD modes identified in the data from December 10, 2007, with all other modes set to zero. These 
modes are the ones subtracted between the two sky images in Fig. [2] The colour scales are in units of Janskys. 



moment i — 27r|u|. At the smallest |u|, a bin width equal 
to that of the primary beam (20 wavelengths) was used. 
At larger |u|, each bin width is increase by 60 per cent, to 
compensate for the decreasing density of visibilities. 



strapped average of the power spectra over all 10 possible 
cross-correlation pairs of the 5 nights. 



It was found that outliers tended to skew the mean 
power in each annulus. Since the median is much more ro- 
bust to such outliers, we calculate the power in each annulus 
(that is, at each angular scale 1) as the median value over 
all frequencies. The error in each annulus is estimated as the 
median of the absolute deviations from the median power, 
weighted by the noise. The final power spectrum is the boot- 



The power can also be expressed in terms of the 
wavenumber perpendicular to the line of sight, k± ~ 
(^/6608) /ilMpc"^, which becomes useful when discussing 
the 3D power. In the 2D case we will continue to use I. 
Since we do not yet include line of sight information, this 
is the power as a function of I with fixed k\\ =0, which we 
denote P(^|A:|| = 0). 
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3 FOREGROUND REMOVAL 

3.1 Singular Value Decomposition 

Foreground removal techniques typically rely on the fact 
that the foreground signal is expected to be much smoother 
in frequency (that is, has fewer degrees of freedom) than the 
reionization signal, which decorrelates on th e order of one 
to a few megahertz (|Bharadwai Alill2005l ). Observations 
of foregrounds around 150 MHz with GMRT have shown 
that the fluctuations in frequency are large enough to make 
polynomial fits ins ufficient to model them (iAli et alTbOOSl : 
iGhosh et al]l2012h . In this work, we instead use a singu- 
lar value decomposition (SVD), which still isolates smooth 
foreground modes but does not make a priori assumptions 
that the foregrounds can be approximated wit h a particular 
functi on. A similar technique h as been used bv I Chang et alj 
and iMasui et al. to clean foregrounds for HI 

intensity mapping at z ^ 0.8, where the relative domi- 
nance of foregrounds over the 21 c m signal is comparable 
to z ~ 8.6. For reionization jLiu Teg mark (2012 ) have de- 
veloped a framework for using SVD modes of a frequency- 
frequency correlation matrix to clean foregrounds at MWA. 

We perform an SVD for each baseline individually on 
the visibilities arranged in a matrix by time and frequency. 
The number of modes is limited by the 64 frequency chan- 
nels. Fig. [1] shows the singular values for the shortest base- 
lines. The spectra of values on a given baseline is generally 
consistent from day to day, but occasionally there are large 
jumps in both amplitude and rate of decline with mode num- 
ber, which are likely due to either RFI or calibration errors. 
In these cases, the noise on the baseline also becomes much 
larger, such that in the final calculation of the power spec- 
trum their contribution is significantly down- weighted. 

A sky image using 8 hours of data from a single night is 
shown in Fig. [2l compared with the same data after the first 
8 SVD modes, shown in Fig.[3l are removed. The overall flux 
is reduced substantially after only a few modes are removed. 
While the sources in the centre of the field are removed 
quite well, the dominate residuals are the point sources near 
the edge of the beam. This is generically true of any fore- 
gr ound subtrac t ion u sed on this data set, as was also seen 
in Paci ga et al ] (|201lh . This is most likely due to beam edge 
effects, the worst residuals being close to the first null where 
the frequency dependence of the beam pattern is most sig- 
nificant. Though there are sophisticated schemes that may 
be able to model point sources while minimizing the impact 
on the 21cm signal (e. g., [Datta et al.l 120101 : IBernardi et al.l 
I2OIII : iTrott et al.ll2012l l. at the angular scales we are inter- 
ested in for this work (i < 2000) the point sources are con- 
fusion limited and contribute in the same way as the diffuse 
background. 

Each night goes through the SVD foreground removal 
separately, and then the cross-correlations are used to arrive 
at a power spectrum using the method described in section[2l 
The spectra for several numbers of SVD modes removed are 
shown in Fig. 21 



3.2 Quantifying Signal Loss 

A general problem with any foreground removal strategy is 
that it is impossible to completely separate the foregrounds 
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Figure 4. Power spectra before and after SVD mode removal. 
The blue line is before any modes are removed. The green, red, 
cyan, and purple lines are for 4, 8, 16, and 32 modes removed, 
respectively. The error bars are from a bootstrap analysis of all 
cross-correlated pairs from the 5 nights of data. Th e solid yellow 
line is the theoretical signal from ljelic et al.l (|2008h . 

from the signal, such that the foreground remo val will likely 
remov e some signal as well. Early work by iNitv anan 
(|2010h used a simple model of an SVD applied to a single vis- 
ibility matrix to show that the signal loss could be calculated 
analytically. Our method of using an SVD for each baseline 
independently is more complex, and we wish to estimate the 
signal loss directly from the data itself. To quantify the sig- 
nal loss, we aim to find the transfer function between the 
observed power Psvd(^) and the real 21 cm power P2icm(^). 
Since the real power is unknown, we use a simulated signal 
as a proxy. This is added to the data before the foreground 
subtraction and the resulting power spectrum after subtrac- 
tion is compared to the input signal. 

The simulated signal we use is a Gaussian random field 
with a power spectrum from CAMI^ scaled to z = 8.6 using 
the linear-regime growth from z = 1.5, and with the ampli- 
tude calibrated to be similar to the expected 21 cm signal 
from EoR. Fig. [5] is an image of the simulated signal as it 
would be seen by GMRT in the absence of any foregrounds or 
noise. The effect of the beam profile on the power spectrum 
is less than 3 per cent for scales in the range 40 < i < 2000, 
and so has a relatively small effect on the result. 

Given a data set x which is the sum of the observed data 
and a simulated signal, the transfer function T{x F{x);i) 
measures how much of the signal survives in the foreground 
filtered data set F(x) as a function of i. While F{x) can 
stand in for any filter method applied to the visibilities x, 
for the SVD we must also specify that the modes removed 
are those calculated from the visibilities x themselves. While 
the transfer function measures the signal loss for a single set 
of data, the power is measured from the cross-correlations 
of those data sets, so the relationship can be written as 

PsVd(^) = T{X ^ F{xyjfP21cm{i). (1) 

Unless it is necessary to be explicit about the mapping T is 
measuring, we will shorten this notation to simply T{i). 

^ Code for Anisotropics in the Microwave Background; 
|http: //camb. inf o| 
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Figure 5. An image of the simulated signal as seen with GMRT 
in the absence of foregrounds or noise, using the same baselines 
and field of view as Fig.s[2]and[3] The colour scale is in units of 
mJy. The solid circle is the FWHM of the primary beam. 
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Figure 6. Transfer function Ti with 4 (green), 8 (red), 16 (cyan), 
and 32 (purple) SVD modes removed, showing the fraction of 
the 21 cm signal that we estimate survives the SVD foreground 
removal. With only 4 modes removed, most of the 21 cm signal is 
expected to survive. However, when 32 modes are removed, about 
20 per cent or less survives, depending on the angular scale. 



There are numerous ways one can estimate this func- 
tion. The most direct way is to cross-correlate F{x) with 
the injected signal, and normalize by the auto-power of that 
same signal. This is written as 

F(data + signal) x signal 



signal X signal 



(2) 



In the ideal case where F{x) removes foregrounds perfectly 
this will equal exactly 1.0. While conceptually simple, and 
used successfully bv lMasui et al.l (|2012l ) for data dit z ^ 0.8, 
we find this estimator of the transfer function to be excep- 
tionally noisy for realistic cases where F(x) leaves residual 
foregrounds. In the case of the SVD, we would expect the 
function to become less noisy as more modes are removed 
and the residual foregrounds decrease, but we are still sig- 
nificantly limited in being able to measure the power. 

An alternative is to subtract the original visibilities, un- 
der the same foreground filter, from the combined real and 
simulated visibilities before cross-correlating with the sim- 
ulated signal. To distinguish it from the previous one, we 
denote this version of the transfer function Ti, which takes 
the form 



[F(data + signal) — F(data)] x signal 
signal X signal 



(3) 



In addition to being much less noisy when residual fore- 
grounds are present, this has the benefit that by subtract- 
ing the original data we remove the possibility of the real 
21 cm signal in the data correlating with the simulated sig- 
nal and biasing the result. If F(x) left the signal untouched, 
this would in principle equal 1.0. However, deviations are 
possible even when F{x) = x. This is due to the fact that 
the cross-correlations with real data in the numerator in- 
troduces RFI masks, noise, and day-to-day variations which 
are not present in the pure signal in the denominator. Thus, 
the transfer function will also correct for these effects, which 
enter at a level of a few per cent. 

We carry out this process of estimating T{i), averag- 
ing over 100 realizations of the simulated signal, after which 
both the mean and the standard deviation are well deter- 
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Figure 7. The Ti corrected power spectra after SVD foreground 
removal. The colours represent 4 (green), 8 (red), 16 (cyan), and 
32 (purple) SVD modes removed. The dotted lines are the un- 
corrected power spectra as in Fig. |4l while the solid lines are the 
power spectra after correcting for the transfer function. For 32 
modes removed, at low i, the corrected power is larger than that 
for only 16 modes removed. As in Fig. 0] the solid yellow line is 
the theoretical signal from ljelic et al.l (|2008|). 



mined, and the error in the mean is small enough that it will 
not contribute significantly to the corrected power spectra 
later. Fig. [6] shows Ti for a selection of SVD filters. While 
the transfer function in principle can depend non-linearly 
on the amplitude of the input signal, we find that the result 
does not change significantly within a factor of 10 of real- 
istic signal temperatures. In the regimes where the transfer 
function does begin to depend on the input temperature, 
the two are anti-correlated; larger signals are more readily 
misidentified as foregrounds by the SVD, leading to a small 
value of T{£). 

The transfer function can be used to determine the best 
number of modes to remove, since as more modes are re- 
moved more of the 21 cm signal will be reduced to a point 
where the additional correction to the signal outweighs the 
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Figure 8. Hermite window in frequency space. Four examples are 
shown with different values of increasing from left to right, with 
arbitrary horizontal offset. Bars indicate the value of the window 
in each frequency bin. In practice should not be smaller than 
the frequency resolution, limiting the /c|| available. The inset plot 
shows the Fourier transform. 



gain from reducing the foregrounds. Fig. [7] shows that cor- 
recting for the transfer function after 32 modes are removed 
gives a weaker hmit on the power than only removing 16 
modes. 



4 3D POWER SPECTRUM 
4.1 Line of sight power 

The power calculated from annuli in visibility space only 
measures the 2D power perpendicular to the line of sight 
(that is, as a function of the multipole moment i or 
wavenumber k±). To find the full 3D power, we must also 
look at the line of sight, or frequency, direction and measure 
power as a function of k\\. While certain forms of foreground 
filters will have a window function that naturally selects a 
k\\, the SVD filter does not have a well defined behaviour 
along the line of sight. The gives us the flexibility of select- 
ing the window function. 

Hermite functions, having the benefit of zero mean and 
a simple Fourier transform, are well suited to select a range 
of k\\. In frequency space, we define a window 



1 - ^ ) exp I 1 



2C 



(4) 



where C is a parameter analogous to the standard deviation 
of a Gaussian distribution, which in this case specifies the 
location of the zeros. This is shown in Fig. [8] for several ^ 
compared to the frequency bin size. This window has the 
Fourier transform 



- exp 



(5) 




We have used the conversion factor r ^ 11.6 /i"^Mpc/MHz 
such that k\\ is in units of hMpc~^ . The normalization has 
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Figure 9. An example power spectrum at = 0.25 MHz (/c|| = 
0A9 hMpc~^) with 16 SVD modes removed (blue line), with the 
effect of each type of transfer function correction for Ti shown. 
The green line corrects for both the Hermite window and the SVD 
subtraction, while the red line reintroduces the Hermite window, 
which agrees quite well with the semi-Hermite correction (cyan). 
The purple line uses the SVD only transfer function. Error bars 
are from include contributions from the transfer function and the 
bootstrap error from the raw power spectra. In this example, the 
three approaches agree quite well, though they can diverge by an 
order of magnitude for other selections of mode subtraction and 
k\\. 



been chosen such that the maximum of h{k\\) is 1, thus pre- 
serving power. This peak in Fourier space, shown in the inset 
of Fig. [H occurs when k\\ = a/2/ (rC) and determines the k\\ 
at which most power survives the Hermite window. Larger ^ 
sample smaller k\\, with the range of possible values limited 
by the frequency resolution and bandwidth. 



4.2 Three approaches to the transfer function 

By applying a Hermite window to the data, we can calcu- 
late the 2D power spectrum at a fixed k\\. There is some 
complication, however, in how we can use a transfer func- 
tion to correct for possible signal loss. The Hermite filter 
by design reduces power on most scales while leaving power 
only at a specific k\\ , and we want the transfer function to 
only adjust for signal lost at the same k\\. Ideally one would 
apply the Hermite filter first to isolate the input power at 
the scales of interest and run the foreground filters on that 
data. If H(x) represents the data set with a Hermite window 
applied, this would measure T{H{x) F{H{x));i). Unfor- 
tunately, the SVD actually depends strongly on information 
in the k\\ direction, which means that F(H(x)) may have a 
much different effect on the power at the chosen length scale 
than F{x). That is to say, the Hermite and SVD operations 
do not commute. 

There are several possible approaches to get around 

this: 

(i) Assume that the transfer function is not strongly de- 
pendent on k\\ , and use T{x F{x);i) from the k\\ = case 
independent of the k\\ selected by the Hermite window. We 
call this the 'SVD only' approach. The k\\ behaviour only 
enters in through calculation of the power spectra after the 
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Hermite window. Any important behaviour of the SVD in 
the A; II direction will not be captured. 

(ii) We can calculate a transfer function for the signal 
loss due to the total effect of both the Hermite window and 
the SVD, T{x H{F{x))]£), and correct for both. We can 
then use an analytical form for the transfer function of the 
Hermite window alone to reintroduce the scale window and 
keep only the power at our selected k\\ . We call this the 'full 
Hermite' approach. 

To find its analytical form, we start with the fact that 
the transfer function associated with the Hermite window 
measures the ratio of the windowed power to the full power. 



^2,, , _ /J^(fc±,fc||)|fe(fc||)Pdfc|| 
^' ~ /P(fcx,fc||)dfc|| • 

If we assume the power spectrum has the form 
1 



P(/C|| 



(6) 



(7) 



both the numerator and denominator of this can be evalu- 
ated analytically. The result is 



7T \ 



(8) 



where erfc[x] = 1 — erf [x] is the complimentary error func- 
tion. Requiring the most steps, this method has more av- 
enues to introduce errors or biases. 

(iii) Apply the Hermite window first to the simulated sig- 
nal. When added to the full data and passed through the 
SVD foreground removal, the larger amplitude of the fore- 
grounds present in the data ensures that the SVD still has 
data at all k\\ to operate on. However, since there is only a 
simulated signal at a specific /cy , the cross correlation with 
the simulated signal when calculating the transfer function 
T (data +i7 (signal) F(data+ if (signal) ); ^) only measures 
the effect of the SVD on that k\\ . We call this the 'semi- 
Hermite' approach. This assumes that the SVD as applied 
to the k\\ limited simulated signal is a suitable proxy for how 
the SVD affects the real signal, given that both the real sig- 
nal and the k\\ limited simulated signal are of significantly 
lower amplitude than the foregrounds. 

Fig. [9] shows a typical power spectrum for a particular 
choice of SVD filter, transfer function, and k\\ , without any 
correction and the resulting spectra after each of the above 
approaches. Differences in each approach illustrate the dif- 
ficulty in finding an unbiased estimator that gives a robust 
result. Fig. [TOl shows the transfer function as a function of 
both k± and k\\ for each of the three possible approaches. 

It can be seen in both the full Hermite and semi-Hermite 
methods that there is a /cy dependence which is not captured 
by the SVD only method, which is constant with k\\ by def- 
inition. It is also notable that the 'full Hermite' approach 
finds Ti further away from 1 than the other approaches for 
zero SVD modes removed. This is likely indicative of a mis- 
match between the amount of power being removed by the 
combination of SVD and Hermite filters and the amount 
modeled by the analytic form. This suggests that in areas of 
(/c^, k\\) space where Ti > 1, this method may overestimate 
the amount of signal present, in turn underestimating the 



21 cm power by failing to fully correct for the signal loss. 
None the less, the full and semi-Hermite approaches agree 
much better with more SVD modes removed. Since the semi- 
Hermite approach seems to capture both the k\\ dependence 
and is relatively well behaved with < Ti < 1, we use it as 
the canonical transfer function. 



4.3 Sampling (k±^k\\) space to get P(k) 

Using the Hermite window to select a fixed k\\ allows us 
to calculate P{l\k\\) and the associated transfer function at 
that k\\. By repeating this for a series of /cy , we can build up 
the full 3D power spectrum. 

Fig. [11] shows the power as a function of both k± and 
k\\ using the semi-Hermite correction, given a series of dif- 
ferent SVD mode subtractions. The power shows a pat- 
tern of lower values towards low k± and high k\\. Fig. [T2l 
shows the same measurements as a function of the 3D 
wavenumber k — \/k±_'^ + k\\^ . Though the SVD is our pri- 
mary mode of foreground removal, the Hermite function 
itself acts as a foreground filter removing the large scale 
structure in frequency space. This is reflected in the points 
where zero SVD modes have been removed. It is clear that 
our ability to remove foregrounds drops off quickly above 
about k ^ 0.5/iMpc~^ Our best limit at 2a is (248 mK)^, 
achieved at {k±,k\\) = (0.11, 0.49) /i Mpc~^, or a total k of 
0.50 /i Mpc~^, with 4 SVD modes removed. At this point, the 
semi-Hermite value of the transfer function was Ti = 0.74, 
meaning an estimated 26 per cent of signal was removed by 
the SVD mode subtraction and Hermite window operating 
on each day in the cross-correlations. If instead 16 modes are 
removed, the limit changes to (319 mK)^ but 55 per cent of 
the signal is lost. Any residual foregrounds, though reduced 
by a much larger fraction than the signal, will also have 
been boosted by this correction, making this measurement 
an upper limit on the actual 21 cm signal. 



5 CONCLUSION 

Using a singular value decomposition as a foreground re- 
moval technique and a simulated signal to quantify the loss 
of a real 21 cm signal the SVD may cause, we have calcu- 
lated an upper limit to the HI power spectrum at z = 8.6 
of (248 mK)^ at A: = 0.50/iMpc~\ The k± component was 
found using the median power in annuli of the {u, v) plane, 
while a Hermite window was used to sample the k\\ direction. 
This is in contrast to our previous work with a piecewise- 
linear filter which operated only in the frequency direction 
and carried with it an implicit k\\ window. 

This limit is dependent on the method one chooses to 
calculate the transfer function between the real 21 cm signal 
and the observed power. Both the k± and k\\ behaviour of 
the foreground filter chosen needs to be taken into account. 
While the semi-Hermite method chosen uses a simulated 
signal with power in a limited k\\ window, and may miss in- 
teractions between the SVD filter and the signal over larger 
k\\ bands, we believe it to give the most reliable estimate of 
the transfer function and a suitably conservative estimate 
on the final upper limit. 

Had we instead used the full Hermite approach de- 
scribed, this limit would have been (260 mK)^. That this 
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Figure 10. The transfer function Ti as a function of both k± and /c|| calculated using the full Hermite (left column), semi-Hermite 
(middle column), and SVD only (right column) methods. From top to bottom, 0, 4, 8, 16, and 32 SVD modes are removed. Slight 
deviations from unity with zero modes removed illustrate the additional effects from RFI masking, noise, and the beam that are captured 
by the transfer function. 



second approach gives a similar value suggests that this limit 
is a fairly robust one. The difference can likely be attributed 
in part to the simplifying assumptions necessary when de- 
riving the analytical Hermite windowing function. We also 
consider the curren t result to be more robust than that re- 
ported previously in lPaciga et al.1 (|201lh . While the previous 
limit was considerably lower, this can be accounted for by 
many factors; the different k scale, the change in foreground 
filter, several minor changes in the analysis pipeline detailed 
in section [21 and most significantly the fact that this is the 
first time a transfer function has been used to correct for 
signal lost in the foreground filter. Without such a correc- 
tion, our best upper limits with the SVD foreground filter 
may have been incorrectly reported as low as (50 mK)^. 

A value of the Hi power spectrum of (248 mK)^ is almost 
an order of magnitude higher than what is generally consid- 



ered physically plausible, and can not constrain reionization 
models with a warm IGM where the spin temperature is 
much greater than the CMB temperature. If we assume a 
cold IGM, a neutral fraction of 0.5, and fully ionized bub- 
bles with uniform radii, the brightness temperature of the 
neutral IGM would have to be at least 540 mK in absorption 
against the CMB to create power this large. 

This limit still compares favourably to others estab- 
lished in the literature which are on the order of several 
Kelvin (e.g., Bebbington 1986; Ali et al. 2008; Parso ns et all 
l201Qh . LOFAR has begun publishing initial results from 
reionization observations, but have so far focused on much 
longer scal es (i ^ 7500) which are not directly compara- 
ble to ours ( Yatawatta et al.ll2013 ). PAPER has new results 
forthcoming (Parsons et al., in preparation) which promise 
to be competitive with those presented here. 
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Figure 11. Power spectra as a function of both k± and k\\ corrected with Ti calculated using the semi-Hermite approach. The topmost 
plot shows the entire k± range without any foreground removal. The four smaller plots show 0, 4, 8, and 16 SVD modes removed on 
the same colour scale for only the lowest few k± bins. The colour scales are in units of log(K^). Compared to the case with modes 
removed, the SVD tends to reduce the overall power by 1-3 orders of magnitude. See also Fig. [TJl which shows the power as a function 
of the total k. 



The SVD procedure could be refined further by a 
baseline-by-baseline accounting of the optimal number of 
modes to subtract or by limiting the field of view on the sky 
to the innermost area of the beam where point source resid- 
uals are minimal, although it is not obvious what effect this 
would have on the signal at small angular scales. Making a 
measurement at larger i would require a more careful treat- 
ment of point sources but is also limited by the fact that 
the SVD is less effective for longer baselines. Regardless of 
the foreground removal technique used, it is likely that ac- 
curately correcting for the any resulting loss of the 21 cm 
signal, and disentangling the 21 cm signal from any residual 
foregrounds, will remain a significant challenge in measuring 
the true EoR power spectrum. 



6 ACKNOWLEDGMENTS 



Special thanks to Eric Switzer and Kiyoshi Masui for many 
useful discussions during the preparation of this work. We 
thank the staff of the GMRT that made these observations 
possible. GMRT is run by the National Centre for Radio As- 
trophysics of the Tata Institute of Fundamental Research. 
The computations were performed on CITA's Sunnyvale 
clusters which are funded by the Canada Foundation for 
Innovation, the Ontario Innovation Trust, and the Ontario 
Research Fund. The work of GP, UP, and KS is supported 
by the Natural Science and Engineering Research Council 
of Canada. KB, JP, and TV are supported by NSF grant 
AST-1009615. 



© 0000 RAS, MNRAS 000, 000-000 



10 G. Paciga et al 



10^ 




Figure 12. Power as a function of the total wavenumber k = -sjhx^ + fcy ^. Each point represents a different (fc_L,/c||) pair; there is no 
binning in k. Colours indicate the number of SVD modes removed; (blue), 4 (green), 8 (red), 16 (cyan), and 32 (purple) are shown. The 
boxed region at 0.5 is shown inset, with nearby points each of the three marked k spread out slightly for clarity. The best limit at 2a 
is (24 8 mK)-^ at 0.50 /iMpc"-*^ achieved with 4 SVD modes removed. The solid line is the predicted 3D power spectrum from ylliev et al.l 
(|2008l ) assuming a 30 mK signal. 
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